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In this paper, a delocalization-localization transition within the superconducting state is explored. 
The symmetries of the Bogoliubov deGennes Hamiltonian endow the two associated superconducting 
phases - the thermal metal and the thermal insulator-, and the critical point between them with 
properties that are in stark contrast with analogous phases in normal systems. Here, for systems 
preserving spin rotational invariance and time reversal symmetry in three dimensions, the transition 
is shown to be of a different universality class from its normal partner by extracting and comparing 
the localization length exponent. The density of states, which may be regarded as the 'order 
parameter' for the field theory describing the superconducting system, is studied for its unusual 
properties about criticality. 



I. INTRODUCTION 

The problem of Anderson localization in normal disor- 
dered electronic systems has been studied for years, and 
has continued to pose intriguing puzzles ^$2\ . Recently, 
the question of whether an analogous transition can occur 
in the behaviour of low-energy quasiparticle excitations 
about a superconducting ground state has come to light 
(see for e.g. Ref. ||.) In fact, the stability of two very 
distinct superconducting phases which may be charac- 
terized by their transport properties - one with extended 
quasiparticle states at the Fermi energy, and the other 
with localized states- has been established Q]. Here, it 
is the transition point between these two phases, and its 
critical properties, that forms the subject of our atten- 
tion. 

While one can afford to ask of the superconducting 
systems the same questions that have been addressed 
in normal systems, the former displays refreshingly new 
physics that shows both conceptual and qualitative dif- 
ferences from the latter. To begin with, quasiparticle 
excitations in superconducting systems do not conserve 
charge, and thus one cannot study the transition through 
charge dynamics. This absence of U{\) symmetry in the 
Bogoliubov deGennes Hamiltonian, which we employ to 
describe the excitations in question, has marked conse- 
quences. As emphasized by Ref. ||, it gives rise to a host 
of new universality classes ,with critical exponents whose 
values are significantly different from those of their nor- 
mal partners. In particular, the density of states(DOS) 
exhibits astonishing features. In this paper, we explore 
these surprises in the context of superconducting systems 
that respect spin-rotational (SU(2)) invariance and time- 
reversal (T) symmetry. 

The issue of a delocalization-localization transition 
within the superconducting state has been a recurrent 
theme (see for e.g. Ref ||). But it is only in recent years 
that analyses of both microscopic and continuum models 
have paid focused attention from a variety of different av- 



enues [HJ4|J^ |Tl[ , and delved into the prospect of making 
this transition realizable in physical systems. To repeat 
the example offered of the dirty d-wave superconductor, 
consider a system of impure superconducting sheets with 
d-wave pairing coupled to one another. At the nodal 
points, one has low-energy quasiparticle excitations, and 
in fact, due to disorder, one even expects a finite DOS 
at the Fermi energy For low interplane coupling 

strength or high impurity concentration, one would ex- 
pect these states to be localized, and upon increasing 
the coupling or lowering the disorder, one could conceive 
of accessing a critical point beyond which these states 
become extended. Here, we study such a transition on 
more generic grounds by exploiting the field-theoretic set- 
up offered by Ref. Q, and through numerical analyses, 
both of which serve to bring out the novel features of 
superconducting systems quite dramatically. 

We review the framework used to describe quasiparti- 
cle excitations in the absence of interactions. We expand 
on some details of the two phases, the 'thermal metal' 
and the 'thermal insulator'. Then, in order to keep our 
paper self-contained, we elaborate on the field-theoretic 
and numerical methods. We proceed to discuss the field- 
theoretic and numerical predictions for the localization 
length exponent v associated with the thermal metal- 
thermal insulator transition, and show that they both 
confirm the existence of a new universality class with 
v < v n , where v n is the corresponding exponent for nor- 
mal systems. We then study the unique properties of the 
DOS at criticality. Finally, we mention the character- 
istics of transitions in superconducting systems besides 
those that respect SU(2) and T symmetry, and possibil- 
ities for experiment. 



II. MODELS AND SYMMETRIES 
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A. The Superconducting Hamiltonian 

Quasiparticle excitations about the superconducting 
ground state are well described within the framework of 
the Bogoliubov deGennes(BdG) Hamiltonian. For the 
spin-singlet paired superconductor, which we focus on 
here, the BdG Hamiltonian has the form 



Ho = Hi + H 2 , 



H 2 



fdfixcUx) ( i ~ iH ^2 Ax)T ~ Ef + V(x)j c a (x\ 



= J d d xd d x'(c\(x)A(x,x')c\(x') + c l (x)A*(x,x')c< [ (x')), 

where c 1 and c are electron creation and annhilation op- 
erators respectively, m the mass , Ep the Fermi Energy 
, V(x) a random potential describing impurities in the 
system, and A(x) a vector potential describing any ex- 
ternal magnetic field. The lattice version of Hq, which is 
more tractable for numerics, has the form 



h.c 



(1) 



From Hermiticity, one requires the condition tj 
while spin rotation invariance requires Ay = Aji. The 
BdG Hamiltonian contains anamolous terms that reflect 
the fact that the excitations do not conserve charge, and 
hence break the associated U(l) symmetry. As a result, 
the BdG Hamiltonian lives in an 'extended particle-hole 
space' with twice the degrees of freedom of normal elec- 
tronic systems. More importantly, the absence of U(l) 
charge conservation endows the BdG Hamiltonian with 
symmetry properties which are completely different from 
normal systems. It is this difference that plays the key 
role in giving rise to properties that are unique to super- 
conducting systems. 

The symmetries of the BdG Hamiltonian have recently 
been studied in the context of mesoscopic systems and 
random matrix theory 0. Within the class of models 
described in Ref ||, the Hamiltonians that we study in 
this paper have SU(2) and T symmetry; they describe 
systems that are singlet-paired and have no spin-orbit 
scattering, thus having spin-rotational invariance, and in 
addition, they have time-reversal invariance. 

Since the quasiparticles conserve both spin and energy, 
one can explicitly write Eq.[l] in terms of conserved quan- 
tities by defining a new set of fermionic d— operators: 



(2) 



Thus one makes a particle-hole transformation on the 
spin-down operators, but leaves the spin-up operators 
unchanged. The Hamiltonian of Eq.^now takes the form 



Hi. »,( $ A ;i U Y<UI, ! d r C-!) 



where time reversal invariance requires fcy and Ay be 
real. Spin rotational invariance now requires 



ayH ij ay = -H* i , 



(4) 



where a y is the standard Pauli matrix. Spin conserva- 
tion along the z-direction is evident from the fact that 
the physical spin is related to the number operator for 
the 'd-particles' by 



Sf = 



{4c 



i 



(5) 



and that H L of Eq.| conserves particle number. 

The Hamiltonian H.l, in principle, may be diagonal- 
ized by solving the eigenvalue equation 



= E 



u{i) 
v(i) 



(6) 



u (j)* 










. -«*(<) . 



where u and v describe the wave-function amplitudes at 
each site. Given Eq.||, one can construct the following 
state 



(7) 



with eigenvalue — E; SU(2) invariance requires a sym- 
metric dispersion about the Fermi energy, and eigenval- 
ues of the BdG Hamiltonian come in pairs (E, —E). 

We now draw attention to the physical situation de- 
scribed by the wavefunctions of Eq.^, and the relation 
they bear to the phase transition of interest. To be- 
gin with, in gapless superconductors, which we focus on 
here, one finds states at and about the Fermi energy. 
These states have a profound impact on thermodynamic 
and transport properties, as discussed previously Hfl|l0|. 
Moreover, the transport properties depend crucially on 
the nature of the eigenstates of Eq]| at E — 0, which 
in turn is determined by the spatial configuration of the 
random potential V(x) and the gap-function A(x). In 
particular, depending on whether these eigenstates are 
extended or localized, one can conceive of two very differ- 
ent superconducting phases - the 'thermal metal' which 
is capable of transporting energy, and the 'thermal in- 
sulator' which cannot conduct energy over large length 
scales. In fact, it has been shown that both phases are 
stable in three dimensions Q]. Thus, one can characterize 
the two phases by the thermal conductivity k. In systems 
with SU(2) symmetry, since the same quasiparticles that 
carry energy also carry spin, the spin conductivity cr s , 
may also be used to describe the two phases, and is re- 
lated to k via an analog of the Weidemann- Franz law in 
the thermal metal: 
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(8) 



C. Hamiltonian for Numerics 



Note however, that as the quasiparticles do not carry well 
defined charge, the Weidemann-Franz law breaks down 
with regards to electrical conductivity. 

The two distinct superconducting phases at hand are 
in analogy with, but quite different from, the metallic and 
insulating phases in normal systems. We study the crit- 
ical properties of the transition between the two phases 
using the formalism and physical set-up described below. 



B. A Field-theoretic Formulation 



The predictive power of the effective field theory of 
the previous section lies in the fact that some results de- 
rived from it are universal to all Hamiltonians satisfying 
the appropriate symmetries. In practice, we find that 
in obtaining the localization length exponent, an analog 
of the tight-binding Anderson model commonly used for 
numerics [ [P3||l7| ] shows crisp data with much less noise 
than other models that we have studied. We focus on 
the same model for the density of states since the locali- 
sation length numerics enables us to identify the critical 
point quite accurately. With reference to Eq.[l], the cou- 
plings take the form 



As with the case of Anderson localization [|19[|2[|, the 
problem of quasiparticle localization within the super- 
conducting state can be described within a field-theoretic 
framework. To briefly review its key features, in pre- 
vious work jj], starting with the Hamiltonian of Eq{j], 
T.Senthil.et.al derived a useful field theoretic action from 
which one can extract a rich variety of properties of 
superconducting systems. Coupling the fermionic de- 
grees of freedom in Eq.[l] to an infinitesimal Zeeman 
field r\ (which acts as a chemical potential for the d- 
quasiparticles) , and employing the replica method to cal- 
culate disorder averages, they obtained an effective ac- 
tion. Fluctuations about the saddle-point of this action 
were captured near two spatial dimensions by a non- 
linear sigma model(7VLtjM) treatment, yielding for the 
final form of the action, 



M 



d d x[— Tr(VU.VW)-r)Tr{U + tf)]. (9) 



A, 




i j , n.n. 

i = j 
otherwis 



i j , n.n. 

i = .] 

otherwise 



(10) 



(11) 



where n.n denotes nearest neighbours, and Va and ViA 
are on-site random variables chosen from a uniform prob- 
ability distribution ranging from —W to +W. We work 
with a three-dimensional cubic lattice described by Eq^, 
and having real couplings with the specific form 



H L = -j=[t®a z 
v2 



(12) 



Here, the cr's denote Pauli matrices, and t and A are ma- 
trices with off-diagonal terms set to unity, and diagonal 
terms taking on values V2Vu and V2ViA- 



Here, 't' is a dimensionless coupling constant that has 
the physical interpretation of inverse spin conductance, 
i.e. j = \cr s . U(x) is a matrix with symplectic Sp(2n) 
group structure, where 'n' is the number of replicas. 

The action Snlo-m given above is quite different from 
the analogous action describing normal systems. The 
field theory , referred to as the 'the principle chiral 
Sp(2n) model', has its first term invariant under the 
global 'rotation' U A^UB, where A\B e Sp{2n) , 
thus possessesing Sp(2n) x Sp(2n) symmetry. The second 
term in the action reduces the symmetry to an Sp(2n) 
symmetry as it only allows invariance under U — > A> UA. 
A knowledge of this symmetry structure proves to be very 
useful in deriving critical properties. 

As done for electrical conductivity in normal systems, 
we employ a scaling theory for the inverse spin conductiv- 
ity t, and analyze the critical point seperating the ther- 
mal metal and the thermal insulator. Specifically, we 
extract the localization length exponent and the unusual 
singular behaviour of the DOS at the Fermi energy. 



D. The Transfer Matrix 

While we obtain the DOS from the Hamiltonian of 
Eq.[l2| by the straightforward process of exact diagonal- 
ization, we extract the localization length exponent by 
means of a transfer matrix formulation that caters to 
Eq.[l2| and tremendously reduces the dimensions of the 
matrices involved in numerical work. 



i+l i 




FIG. 1. A quasi 1-dimensional set-up for transfer matrices. 
In-going and out-going wavefunction amplitudes at either end 
are denoted by ipi„ and ipout respectively. 
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To describe its principle, consider a quasi 1- 
dimensional strip in d-dimensions with cross-sectional 
area L d ~ x , and in-going and out-going states at either 
end of this long strip as shown in Fig.|l|. One might for- 
mally obtain the scattering matrix S for the in-going and 
out-going states using the definition 



S = lim exp ■ 



iHT 



(13) 



The scattering matrix, in the specific basis of in-going 
and out-going states can be written in terms of reflection 
and transmission matrices V and 't' respectively fig]: 



S = 



r t' 
t r> 



(14) 



Given this form of the scattering matrix, one can easily 
derive the transfer matrix T which has the property 



ip L = Tip R . 



(15) 



The symmetries of the BdG Hamiltonian imply that the 
transfer matrix too has very specific symmetry proper- 
ties which distinguish it from those of normal systems. 
Of late, these symmetries have been explored on group 
theoretic grounds [p|. 

The transfer matrix T can be constructed by multiply- 
ing a set of transfer matrices, each appropriate for a slice 
of the strip shown in Fig.0: 



-1,R 



N 



T = JjTi. 



(16) 



(17) 



The form of the transfer matrix that we use for numer- 
ical calculation does not make the symmetry of the BdG 
Hamiltonian manifest, but it is tailored specifically for a 
tight-binding Hamiltonian such as the one described in 
Eq.||. We begin by writing the Schrodinger EqJq as a 
difference equation 



AiDi + B 



Bii-iDi-i — EDi, 



(18) 



where Di denotes the wave-function amplitudes on each 
slice in the spinful d-quasiparticle eigenbasis, and A t and 
Bi i+i are 2L 2 x 2L 2 matrices of the form 



A; = 



B 



(19) 



(20) 



where the tij and A^'s are now matrices coupling slices 
T and 'j' in the manner described by Eq.^. To obtain 
the transfer matrix, we rewrite Eq.ftS as 



A+i - B~t +l (E - Ai)Di - /;,,';/-•, ,. (21) 

For our model, the inter-slice coupling has the simple 
form 



Bi. l+1 = A=[I® <J X +I®a z ], 
v2 

which satisfies the special property 



Bi,i+1 — B. 



l 

i,i+V 



(22) 



(23) 



Using Eq.|2l], we are now in a position to define a transfer 
matrix as follows: 



A+i 
Di 



T, 



where T has the relatively simple form 



Ti 



Bi : i+i(E 
I 



A,. 







(24) 



(25) 



and the multiplicative property of EqjlT]. With this 
transfer matrix at hand, which very closely resembles 
standard ones used for the Anderson model |l^Jl^], we 
can extract the localization length for different values of 
energy E and disorder W in a manner completely analo- 
gous to the numerical treatment of the Anderson model. 
As we are interested in the behaviour of states at the 
Fermi energy, we set the energy E, to zero. 

The procedure for extracting the localization length 
is quite standard, and has been elaborated on in great 
depth in many works 0,0. But to briefly outline the 
method, one begins with an orthonormal basis of vectors 
di (0) in the space of the transfer matrix, representing the 
right-most states in FigJIJ One then assumes that for a 
given disorder strength W, and width L, there exists a 
set of eigenmodes Wi(L, W) that describes typical eigen- 
modes for a quasi 1-dimensional system of N slices, where 
'N' is large enough to represent the average behaviour of 
the random disorder, and that these modes decay or grow 
in magnitude as exp(±7iiV) upon multiplication with the 
corresponding transfer matrix T = ]X =1 T. One further 
assumes that when the set of basis vectors di(0) , each of 
which may be represented as a linear combination of the 
Wi(L,W), is multiplied by T, the resulting vectors ^(1) 
(where the dependence on L and W is implicit), are each 
composed of the appropriate linear combination of modes 
Wi(L, W) now weighted by the corresponding growth fac- 
tors exp(±7;7V). The length 'N' is numerically restricted 
by the exponentially growing magnitude of the vectors v. 

Since all vectors have now grown fastest along the 
fastest growing mode, say w\{L 1 W), the magnitude and 
direction of any one of these modes, say Ui(l), are 
exp7iiV and along W\(L, W) respectively. Projecting 
out the component u!(l) from the next vector ^(l) 
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gives a resultant vector 02(1) whose magnitude is roughly 
exp7 2 -/V. Thus, by such a 'Gram-Schmidt' orthogonal- 
ization procedure for the whole set of vectors vl(l), one 
obtains a set of orthogonal vectors 61(1) with associated 
Lyapunov exponents 



7,(1) = 



ln|o1(l)| 
N 



(26) 



which give the characteristic inverse localization lengths 
associated with each mode. To reduce computational 
effort, we consider just the positive Lyapunov exponents 
corresponding to exponentially growing states. As we are 
restricted in our length size N, to obtain a fair estimate of 
the typical 7i's, we repeat the procedure of transfer ma- 
trix multiplication, now using as our initial basis vectors 
the normalized set Oi(l) which more or less point along 
the ideal basis Wi(L, W). An average value 2y=i 7»(i) 
obtained from M such iterations provides the desired es- 
timate of the ideal 7i's. We associate the characteris- 
tic localization length X(L, W) with the slowest decaying 
mode, and thus with the inverse of the smallest positive 
Lyapunov exponent j m in(L, W). 

In the quasi 1-dimensional case, all modes are ex- 
ponentially decaying or growing (corresponding to in- 
going or out-going states respectively) since even the 
slightest disorder is enough to localize states. But in 
the 3-dimensional limit, where the 2-dimensional cross- 
sectional area becomes large, we know that the modes 
ought to experience a transition from extended to local- 
ized behaviour as a function of disorder. To determine 
the critical disorder strength W c for this transition, the 
nature of the modes in 3-dimensions and the localiza- 
tion length for an infinite size system of given disorder 
strength, one can use a finite size scaling analysis of the 
quasi 1-dimensional system. The scaling function that 
we will use to do so is the dimensionless parameter 



k(L,W) 



X(L,W) 
L ' 



(27) 



We now turn to the critical properties of the phase tran- 
sition between the thermal insulator and thermal metal. 



III. CRITICAL PROPERTIES 

As in disordered electronic systems, we have seen in 
Scc. [II B| that the thermal metal, the thermal insulator 
,and the critical point seprating them may be charac- 
terized by their transport properties. The replica field 
theory of Eq|| with its dimensionless coupling t, pro- 
vides a powerful means of studying this transition. An 
analysis of the action in Eq.|| shows that in 2 + e di- 
mensions, where e = 1 for our system, an unstable fixed 
point t c describes the critical point between the thermal 
metal and thermal insulator. One can study the scaling 



behaviour of t with system size L(see Fig.p]) by deriving 



a form for the scaling function [3(t) — d ^ L , which we 
present explicitly in the next section. Near i c , for t > t c , 
the coupling t grows larger with L, and thus exhibits a 
stable thermal insulator, while for t < t c , a smaller and 
smaller value of t with increasing length scale signals a 
thermal metal. 



THERMAL METAL 



THERMAL INSULATOR 







t 



FIG. 2. Behaviour of coupling constant 't' with increasing 
system size,L. 



A. Localization Length Exponent 



The (3 function for the Sp(2n) sigma model of the ac- 
tion of Eq.|| can be found in Ref. , and it is given to 
cubic order('two' loop) in coupling t by 



0(t;Sp(2n)) 



-et + (2n + 1) + i(2n + l) 2 t 3 + 0(< 4 ), 



(28) 



where e = d — 2, and n denotes the number of replicas. 
In the limit n — > 0, at the critical point where the (3 func- 
tion vanishes, its derivative gives the inverse localization 
length exponent: 



0(e 3 ). 



(29) 



In contrast, normal systems with time reversal symme- 
try and spin-rotational invariance may be identified with 
the Sp(4n)/Sp(2n) x Sp(2n) model of Ref. @, and its 
associated f3 n function has the form 

f3 n = -et + (An + l)t 2 + (8n 2 + 2n)t 3 + 0(t A ). (30) 



The localization length exponent derived from Eq.[ 
the value 



1 



0(e% 



has 



(31) 



quite different from the value of v in Eq.g9|. The numer- 
ical evidence to follow supports the field theoretic result 
that the value of the localization length exponent v for 
the superconducting system with T and SU(2) is consid- 
erably lower than the analogous exponent v n for normal 
systems. 
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1. Numerical Treatment 

The standard numerical technique that we employ for 
extracting the localization length exponent i>, shows that 
in three dimensions it takes on the value 1.15 ± 0.15. 

The hnite-size scaling technique can be summarized as 
follows: scaling arguments require that the only relevant 
length scale in the system be the localization length £,(W) 
of the infinite sized system, and thus we have 



A(L, W) 



= A(L,W) = h( 



aw) 



(32) 



where h is a scaling function yet to be determined. Close 
to the critical point W c , we have the localization length 
£ of the infinite system behaving as 



£ ~ \W - W c 



(33) 



which means that the argument 'a;' of h(x) in Eq.|32| blows 
up at the critical point. However, A is well-behaved and 
finite; this is only possible if we have 



lim h(x) = const., 



(34) 



where the constant refers to independence with respect 
to L for large L. Thus, the critical value A c is common to 
all sufficiently large system sizes. Using Eq.^3[ we rewrite 
Eq.B2 as 



In A(L, W) = f[L-(W-W c 



(35) 



Linearizing the function 'f ' about the critical fixed point 
(W c ,lnA c ) yields 



lnA(L,W) = lnA c + A(W - W C )L^ U . 



(36) 



To procure the value of v, we use an iterative procedure 
|l7f which is equivalent to the widely used procedure of 
Rcf. (f| of performing a least square fit to obtain actual 
values of £,(W). We rewrite Eqj3^ as 

In A = AL 1/v W - (- In A c + AL 1/v W c ) 
= a(L)W-b(L,A c ,W c ), 

assume an initial value for the critical point (W c ,hiA c ) 
from Fig.^, obtain the functions a(L) and b(L,A c ) by 
curve-fitting, extract A and v, determine the value of the 
critical point thus obtained, and repeat the procedure till 
convergence is achieved (which happens rather quickly). If 
scaling is valid, we can collapse our data onto curves of 
mA(W,L) vs L X I V {W — W c ). 

In our transfer matrix calculations, we choose sys- 
tems whose cross-sectional areas have linear dimensions 
of L = 4, 6, 8, 10, and our transfer matrices have dimen- 
sions AL 2 x AL 2 with the given values of L. We choose the 
number of transfer matrices to be multiplied together by 



ensuring that the basis vectors do not grow upto a mag- 
nitude greater than 10 5 upon being multiplied by the set 
of transfer matrices. We utilize a total of 2000 slices in 
the quasi 1-dimensional system for each value of L and 
W. 



L = 4 *; 6 +; 8 o; 10 X. 



^ t ~ - _ 

" St--.- -£---«- 



3.5 4 4.5 5 

disorder.W 



FIG. 3. The scaling function A(L, W) plotted as a function 
of disorder for different system sizes. 



Fig.|| shows the plots for A(L, W) as a function of dis- 
order for L = 4,6,8,10. For fixed disorder W, an in- 
creasing A(L, W) with increasing system size L indicates 
the extended regime, while a decreasing A(L, W) shows 
that the system is in the localized regime. In comparison 
to the Anderson model for normal systems, our simula- 
tions require a much smaller number of transfer matrices 
for relatively noise-free data, and we believe that this re- 
ally is a consequence of the relatively low critical disorder 
strength. 
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FIG. 4. The scaling function In A(L, W) plotted as a func- 
tion of disorder for different system sizes. 



Fig.y shows the data for the iterative procedure which 
jives the value for the localization length exponent 



v= 1.15±0.15. 



(W , In A )= (4.18, -0.46); v = 1.13; L= 4 *; 6 +; 8 o; 10 X 



(37) 



L 1,v (W-W c ) 



10 12 



FIG. 5. Collapse of data for lnA(L, W) given the value 
v — 1.13 for the localization length exponent 



Finally, Fig.|| indeed demonstrates clean data collapse 
close to the critical point. 



worth employing them to study this novel phase transi- 
tion and analogous ones in superconducting systems with 
other symmetries. 



B. Density of States 

The quasiparticle DOS in dirty superconducting sys- 
tems exhibits some of the most stunning differences from 
normal systems. In normal systems, both in the Ander- 
son metal and the Anderson insulator, i.e., in the absence 
of interactions, the DOS remains a smooth continuous 
function across the Fermi energy. In contrast, in gap- 
less superconductors that respect SU(2), well within the 
thermal metal, quantum interference effects cause a sin- 
gularity at the Fermi energy which manifests itself as a 
\f~E cusp in 3-dimensional systems JTo[ . Deep in the ther- 
mal insulator, the density of states exhibits a power-law 
that vanishes at the Fermi energy with the form p ~ 
where a = 1 for systems possessing time-reversal invari- 
ance ||[l0]]. About the critical point, the DOS once again 
shows power-law singularities, the details of which we 
discuss below. The curious form of the DOS (shown in 
Fig.^) has profound impact on thermodynamic proper- 
ties, and in particular, manifests itself in quantities such 
as specific heat and spin susceptibility. 





a. w « w. 



b. w = w 



c. w » w. 



FIG. 6. Behaviour of DOS a) Well within the thermal 
metal b) At criticality c) Well within the thermal-insulator. 



2. Summary of Results 



Both field theory and numerics concur with the fact 
that the localization length exponent in the supercon- 
ducting systems is significantly lower than that of their 
normal partners, clearly indicating a new universality 
class. An e expansion in d = 2 + e dimensions of the 
action in Eq|| shows that the superconducting system 
has a localization length exponent v of (e + y + 0(e 3 )) -1 
in contrast to a v n of (e + 0(e 3 )) -1 for normal systems. 
The field theory would thus predict v = | versus v n = 1 
in 3-dimcnsions. In comparison, one obtains the numeri- 
cal estimate v = 1.15 ± 0.15 versus v n = 1.54 ± 0.08 Jl|] 
in 3-dimensions. 

We must remark that the system sizes and computing 
power utilized in our numerical studies were relatively 
low compared to the current cutting edge procedures. 
As a lot of work has gone into refining techniques with 
regard to normal systems(for e.g. Ref pi), it is well 



1. Discussion of Critical Behaviour 

The field theoretic action of EqJ| not only offers con- 
crete predictions for the DOS, if one were to use Weg- 
ner's analogy with magnetic systems p9| , it provides an 
intuitive picture for the behaviour about criticality. To 
elaborate, the quasiparticle DOS at the Fermi energy, 
which also gives a measure of the magnetization, acts as 
the 'order parameter' of the field theory ||. It is given 
by 

p = lim ^ < Tr(U j + U) >, (38) 
n— >0 An 

where po is the bare DOS, and n the number of replicas. 
The field rj, which has units of energy E, couples to the 
DOS in the action of Eqj9| and might be equated with 
the magnetic field in the magnetic analog. 
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With a little indulgence, one can go further with par- 
allels between the field theory and the magnetic systems, 
as first suggested for normal systems [[L9) : 



F.T. of SnlctM 


Magnetic Systems 


Distance from criticality 

A _ w-w c 

Wr. 


Reduced Temperature 
t 


DOS, p 
Energy, E 

v — ie. 

dR 


Magnetization, m 
Magnetic Field, h 
Magnetic Susceptibility,x 



In normal systems, the analogy is clouded by the fact 
that the DOS is a continuous function of energy and dis- 
order respectively. One can reconcile with this if the DOS 
obeys a power-law form with exponent zero, and in fact, 
one can show this to be the case on field theoretic grounds 
Ull. But in superconducting systems, as we shall see, the 
analogy goes through in quite a striking manner, with a 
whole slew of nontrivial critical exponents: 



p(A,E = 0) 




p(A = 0,E) 


~|25|* 


xHA| 


-7 


£~|A| 


-v 



where we have used the notation of the tables 
above, and £ is the localization length and v the asso- 
ciated exponent described in the previous section. 

In order to derive expressions for (3, 6 and v, we start 
with the 'free energy density' f, obtained from the action 
of Eq| : 



/ 



1 lnZ n 
— -7 hm 

L d n-> o n 



dUe 



(39) 



(40) 



where dU , the integral volume element, takes into account 
the symplectic group structure of the matrices in the ac- 
tion SnloM, n denotes the number of replicas, Z denotes 
the partition function of single system, and the overbar 
above Z n refers to the average over disorder. Near criti- 
cality, f s , the singular part of the free energy density, is 
expected to scale as follows: 



f s (A,E)=C d f(ti y E), 



(41) 



where 'y' describes the scaling form of E, and Eq.|33| 
gives the behaviour of the correlation length £ ~ |A| _l/ . 
Differentiating the free energy with respect to E results 
in the following form for the DOS: 



p(A,E)=C d+y F ± 



(42) 



where F + corresponds to behaviour for A > 0, and F 
for A < 0. To obtain /3, we set E = 0, and compare the 



form of the resulting order parameter p(A, E 
above table, yielding 

(3 = v(d-y). 



0) in the 
(43) 

Here, we require F~(Q) to be finite, and F + (0) = 0. 

To extract S, we impose the physical constraint that 
p be well-behaved and finite at criticality. This requires 
that p{A — > 0, E) be independant of the diverging corre- 
lation length, and thus yields 



o y 

Finally, taking a derivative of p in Eqj4^ 
E gives us the following expression for 7 

7 = v(2y - d). 



(44) 

with respect to 
(45) 



2. Results 

To obtain estimates of critical exponents from field 
theoretic results in 2 + e dimensions, we use the value 
(e+ 3") _1 +0(e 3 ) obtained for v in the previous section, 
and the value \ = § + C(e 3 ) from Ref. §]. Eq|f|4| 
then enable us to determine the critical exponents (3 and 
7 via the relationship Specifically, in the case of 

3-dimensions, substituting the value e = 1, we obtain the 
rough estimates y — 12/5, 6 = 4, v = 2/3, (3 — 2/5 and 
7 = 6/5. 

Shifting our focus to numerical results, the method of 
exact diagonalization reveals that the superconducting 
system at hand does indeed show the novel singular be- 
haviour in the DOS at the Fermi energy. In the data 
shown below, we have once more modelled the supercon- 
ducting Hamiltonian after Eq|^ using periodic boundary 
conditions. Systems of linear dimension 'L' have required 
matrices of dimension 2L 3 X2L 3 , and we have explored 
system sizes with linear dimensions L = 4,6, 8. 





0.12r 
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FIG. 7. Progression of DOS as a function of energy for 
various disorder strengths. 



FIG. 9. DOS at E=0 as a function of disorder strength. 



Fig.^ shows the progression of the behaviour of the 
DOS with increasing disorder strength. As seen in the 
last panel of Fig.Q, the DOS shows a power-law behaviour 
of the form p ~ \E\ about the Fermi energy , Ef, consis- 
tent with expectations for the thermal insulator. 



U 6; DOS data from 1000 matrices; W - 4.24 



\ / 
V / 
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\ I 
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\ i 
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i ' 
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Energy 



FIG. 8. Zoom of DOS about E = for W = 4.24. 



Fig.|| shows a zoom of the DOS about Ep for disorder 
close to the critical strength Wc > for which we have an es- 
timate from the localization length study of the previous 
section. One can easily discern that the DOS plummets 
down quite markedly, and docs indeed exhibit singular 
power-law behaviour. 

L= 4; DOS data from 1728 matrices 



3 4 5 

Disorder, W 



A plot of the DOS at the Fermi energy(Fig^) shows 
that even relatively small system sizes provide numerical 
confirmation of the fact that p(E = 0) acts as the order 
parameter for the field theory of Eq.|| the DOS at Ep is 
finite for low disorder, and it slowly drops to zero beyond 
a critical disorder strength. As discussed in the previous 
section, one would in fact expect the DOS for an infinite 
sized system to behave as p(A, E = 0) ~ |A|", where A 
is the distance from criticality within the thermal metal. 
Scaling arguments for extracting (3 require 



PL (A,E = 0) = \AfY(LA»), A<0, 



(46) 



where pz, is the DOS associated with a system of linear 
dimension L, and Y is a scaling function. One can rewrite 
the above equation in a form more conducive to numerics 
as follows: 



p L (A,E = 0)=L-»Y(AL»' 



(47) 



where Y is yet another scaling function with limiting be- 
haviour Y(x — > oo) = \x\P, reproducing the required de- 
pendence of p(A, E — 0) on A for infinite system size. 



Data for p(E=0,W) and L = 4 *; 6 +; 8 
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FIG. 10. DOS at E=0 for various disorder strengths W, 
and system sizes L 



Fig.|Ty exhibits the plots of pi as a function of disor- 
der for different system sizes L, and we make use of this 
data to procure the value of in Eq.^7]. To extract j3, 
wc perform a fit taking W c , v and (3 as variable param- 
eters. Exploiting the universal nature of the function Y 
in Eq.^, we find the appropriate values of p obtained 
by linear interpolation for a given set of system sizes and 
fixed argument in Y, and plot these on a log-log scale 
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versus system size; the slope for a linear fit of such a set 
of points determines (3. The actual value of (3 is obtained 
by performing the above procedure for different values of 
the argument of Y and taking the average of the j3's thus 
obtained. 



W =4.67; v = 1 .25; ft = 0.1 5; L = 4*; 6+; 8 o 



% * 



-2 2 

L 1 '"A 



FIG. 11. Data collapse of p(E = 0,W) given the form 
Eq.p7|, and the critical exponent /3 = 0.15. 



The set of values = 0.15, W c = 4.67 and v = 1.25 re- 
sult in the data collapse shown in Fig.|ll|. In comparison, 
as remarked at the beginning of this section, the field the- 
oretic result predicts that j3 = 0.4. Once more, as in the 
case of the localization length exponent, we comment on 
numerical accuracy; other simulations using exact diago- 
nalization, for instance, those catering to specific physical 
situations [pof , have used larger sytems sizes and number 
of realizations which would be well-worth employed here. 
However, the above numerics conveys quite clearly that 
the DOS at criticality exhibits a power law suppression 
about the Fermi energy, and that p(E = 0, W) acts as an 
'order parameter' with a non-trivial exponent (3 in sur- 
prising contrast to normal systems which have p(E, W c ) 
smooth about Ep, and a vanishing exponent (3. 



insulator are both stable only in 3-dimensions, and the 
transition cannot occur in lower dimensions where quasi- 
particle excitations are always found to be localized at 
large enough length scales. Due to the absence of time- 
reversal symmetry, the Hamiltonian of Eq.[l] now has 
imaginary couplings. As described in detail in previous 
work pc| ], in this symmetry category, the pinned vortex 
state of a type II s-wave superconductor appears to be a 
fine candidate for exhibiting the thermal metal-thermal 
insulator transition. Low energy quasiparticlcs bound 
to the core of the vortices can tunnel from one vortex 
to another, and as the magnetic field is increased, the 
density and tunneling strength also increase. It is con- 
jectured that there could exist a critical magnetic field 
i/ c 4, within the vortex phase at which the low energy 
states can permeate through the medium to form ex- 
tended states. Field theoretic methods have analyzed 
the properties of these systems, and in particular, have 
shown that in parallel to the case with T, the density of 
states at the critical point has the power-law behaviour 
p ~ |-E| C / 2 , where E is the energy and e = d — 2. 

The presence of spin-orbit scattering or of triplet pair- 
ing introduces new ingredients. It breaks spin-rotational 
symmetry, and as in the case of normal systems, field the- 
oretic arguments predict the presence of a delocalization- 
localization transition not only in 3-dimensions but also 
in 2-dimensions pj.In Ref. [jO], T.Senthil et.al. cast 
their Hamiltonian in terms of Majorana fermions, and 
their formulation, among other things, is highly con- 
ducive to numerical work. They discuss possibilities for 
experiment, and mention heavy fermion systems where 
spin-orbit scattering is prominent, as a physical realiza- 
tion. 

Finally, the tantalizing prospect of a Hall effect in 
superconducting systems has been explored in systems 
with other symmetries as well, i.e., those with SU(2) M. 
Ref H shows that a superconductor with unconventional 
d x 2_ y 2 + id xy pairing symmetry is capable of exhibiting 
a phase with non-vanishing spin and thermal Hall con- 
ductances. Indeed, as in normal systems, sophisticated 
methods such as employing the network model and super- 
symmetric spin chains, have shed light on these systems 
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IV. OTHER SYSTEMS 



V. EXPERIMENT 



We have studied the thermal metal-thermal insulator 
transition in superconducting systems with SU (2) and 
T, and discussed the dirty d-wave superconductor as a 
possible physical realization. Superconducting systems 
with other symmetries too promise such a transition. 

The properties and phases of superconducting systems 
with spin-rotational invariance, but no time reversal in- 
variance are rather similar to our case which preserves 
both symmetries. The thermal metal and the thermal 



The transport properties of normal systems have been 
probed in great depth, and now we see that supercon- 
ducting systems could potentially offer an equally rich 
range of experiments in the thermal metal-thermal insu- 
lator transition. In the previous section, we have men- 
tioned a variety of experimental candidates for study, 
such as the dirty d-wave, the type II s-wave, heavy 
fermion and other superconducting systems. These sys- 
tems must share the feature of gapless superconductivity; 
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one requires states at the Fermi energy since it is these 
states that determine transport properties, and distin- 
guish the thermal insulator from the thermal metal. As- 
sociated with each system, a tunable parameter such as 
disorder or magnetic field ought to be able to access the 
phases. As was previously discussed in detail |T(J, the 
type II s-wave superconductor in the vortex phase offers 
promise as a likely candidate for observing the novel tran- 
sition since in principle, one need only tune the magnetic 
field, and generally the vortex phase exists over a large 
range of field. 

As seen earlier, in contrast to normal systems, the den- 
sity of states shows singular behaviour about the Fermi 
energy for both phases and at the critical point. We saw 
that for systems with T and SU(2), it obeys a power-law 
behaviour of the form p(E) ~ \E\ a , where a — \ well 
within the thermal metal, a = 1 deep in the thermal 
insulator, and field theory predicts a = | at the transi- 
tion, with e — I for 3-dimensional systems. This singular 
behaviour ought to be reflected in thermodynamic quan- 
tities such as specific heat and spin susceptibility, and 
in tunneling experiments. In particular, the temperature 
dependence of the thermodynamic quantities would have 
a form C ~ T 1+a for the specific heat, and x ~ T a for the 
spin susceptibility. Recent experiments of cuprate super- 
conductors have observed a suppression of the specific 
heat close to the Fermi energy pjj . However, measur- 
ing the differing behaviours to detemine the phase might 
prove tricky, especially since we have neglected various 
effects such as interactions which could come into play. 

For systems that preserve spin-rotational invariance, 
quasiparticle excitations about the superconducting 
ground state not only conserve energy, but also spin; 
the spin-conductance <j s could be employed to determine 
whether the system inhabits the thermal metal or ther- 
mal insulator. In the thermal metallic phase, a mag- 
netic field gradient would cause the spinful quasiparticle 
to diffuse across the sample, while in the thermal insula- 
tor, they would be unable to conduct spin. A variety of 
refined spin-injection techniques have been developed in 
semi-conductors to measure spin dynamics (see for e.g. 
Ref. [^2)), but by no means would it be a simple task to 
cater these experiments to superconductors. 

We believe that by far, thermal conductivity measure- 
ments would offer most promise in probing the thermal 
metal-thermal insulator transition. For all superconduct- 
ing systems with their differing symmetries, the thermal 
conductivity k, distinguishes the two phases in that the 
ratio k/T tends to a finite constant in the thermal metal, 
and to zero as T — ► in the thermal insulator. Along the 
lines of previous experiments p3[ , it would be extremely 
interesting to observe the transition by applying a small 
thermal gradient across a superconducting sample and 



measuring k/T as a function of a tunable parameter. 
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